Undulation instability of epithelial tissues 
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Treating the epithehum as an incompressible fluid adjacent to a viscoelastic stroma, we find a 
novel hydrodynamic instability that leads to the formation of protrusions of the epithelium into 
the stroma. This instability is a candidate for epithelial fingering observed m vtvo. It occurs for 
sufficiently large viscosity, cell-division rate and thickness of the dividing region in the epithelium. 
Our work provides physical insight into a potential mechanism by which interfaces between epithelia 
and stromas undulate, and potentially by which tissue dysplasia leads to cancerous invasion. 

PACS numbers: 87.19.R-,47.20.Gv,87.19.xj 



Interfaces between epithelial tissues and stromas often 
present different degrees of undulations. In pre-cancerous 
abnormalities of epithelial tissues — called dysplasia — 
such undulations are often especially pronounced and can 
evolve into long fingers that extend into the stroma IJ . In 
a stratified epithelium, an important indicator of tissue 
dysplasia is the thickness of the region in which cells di- 
vide. While in healthy epithelia only the cells directly at 
the basement membrane divide, cell division in dysplas- 
tic tissues takes place in a larger domain, and in severe 
cases throughout the entire epithelium. 

The instability of monolayered epithelia has been mod- 
eled as the result of a buckling phenomenon [2]. Other 
studies have used the framework of nonlinear elasticity to 
describe the instabilities in growing tissues [3] . As moti- 
vated in earlier work yij and shown experimentally [H [6] , 
tissues behave as viscous liquids on long time scales. This 
is illustrated for example by the existence of surface ten- 
sion at tissue boundaries [71110]. Theoretically, viscous 
descriptions have already been applied in other contexts 
of tissue growth [11] . Here, we propose that the fingering 
of a stratified epithelium originates from viscous friction 
effects driven by cell division. We treat the epithelium 
as a viscous fluid lying on top of a viscoelastic stroma 
(Fig. [iJ. As the epithelium consists mostly of cells, the 
stroma is made of a network of collagen and elastin fibers, 
constantly remodeled by fibroblasts present at low den- 
sities [12]. On short time scales, this network responds 
elastically to deformations, but its constant remodeling 
by fibroblasts allows the tissue to fiow on long time scales. 
A qualitative understanding of the full viscoelastic pic- 
ture can be gained by interpolating between the results of 
the elastic and viscous regimes. In this letter, we present 
these two limit cases. 

For the sake of simplicity, the epithelium and the 
stroma are each considered incompressible. In this case, 
the continuity equation for the epithelium reads daVa — 
kp, where kp is the global production rate of cells, taking 
into account cell division and apoptosis. The associated 
constitutive relation is that of an incompressible fluid 



with shear viscosity rj [TB]: a'^p — rj^daVp -I- dpVa). Here 
the total stress tensor aa/s has been split into a dynamic 
part a'^p and a velocity-independent part —pcSap, where 
Pc is the tissue pressure. The system of equations de- 
scribing the epithelium is completed by the force-balance 
condition daCTap = 0, which leads to: 

rjdadaVp -I- rjdfikp - dpPc = 0. (1) 

Similarly, for an elastic stroma, we obtain: 

^ dadaup - dpps = 0, (2) 

together with daU^ = 0, where is the displacement 
field, /i the shear modulus and Ps the pressure. 

Boundary conditions are as follows. The stress van- 
ishes at the apical surface of the epithelium, taking into 
account the Laplace pressure due to the epithelium apical 
surface tension 7a. At the bottom of the stroma, the dis- 
placement vanishes. At the epithelium-stroma interface, 
the normal component of the velocity is continuous and 
the normal component of the displacement of the stroma 
is equal to the variation of the interface location. The dis- 
continuity of the normal component of the stress is given 




FIG. 1: Schematic representation of a stratified epithelium 
sitting on top of its underlying stroma. The arrows represent 
the qualitative profile of the cell-velocity field driven by cell 
division. The epithelium extends a fingering protrusion into 
the stroma, driven by viscous shear stress. 
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by Laplace's law with interfacial tension 71. Finally, the 
tangential components of the stress are continuous and 
equal to a finite surface-friction term with coefficient 

The physical origin of the instability discussed in this 
work can be qualitatively understood as follows. Con- 
sider a fingering protrusion of the epithelium into the 
stroma and assume for simplicity that cell division occurs 
over the entire height of the protrusion (Fig. [T]). The di- 
viding cells create a flow in the epithelium. Since the cells 
above the finger have more dividing layers underneath 
them than their neighbors, they flow toward the apical 
surface faster than the cells in the adjacent regions. This 
results in a shear flow of cells within the epithelium. The 
associated shear stress builds up pressure at the bottom 
of the finger, favoring the development of the protrusion. 

Let us now discuss the solution of our model for the 
flat, unperturbed epithelium-stroma interface. Here and 
in the following, we make the assumption that, due 
to the lack of nutrients and growth factors away from 
the stroma, the overall cell production decreases ex- 
ponentially over a length scale / with increasing dis- 
tance Az from the epithelium-stroma interface: fcp — 
k exp{—Az/l) — kg [T7]. When the interface is flat, the 
cell velocity and pressure in the epithelium read: 



/cm 1 — exp — 



z- L 



-ko(z-L), (3) 



= 2?7 ( fcexp ( -^j^ ] - ko 



(4) 



where the origin of the z coordinate is at the bottom 
of the stroma, and L is the stroma thickness. The 
height of the epithelium H is determined from the con- 
dition that the cell velocity vanishes at its apical sur- 
face. Together with Eq. ([3|, this condition reads ko — 
kl/H{l — exp(— The deformation of the stroma 
vanishes everywhere (m° = 0). 

We now address the question of the stability of the sys- 
tem under a small perturbation. Since we do not expect 
the origin of the instability to depend on dimensionality, 
we consider the case of a system translationally invariant 
in the y-direction, with a perturbation of the epithelium- 
stroma interface of the form 5h{x, t) = SHq exp(wt -l-igx). 
In the linearized system of equations, the solutions for 
the perturbations all take this form. Eq. ([T]) then reads: 



0, 



(5) 



together with the continuity equation daSva = kp{z — 
Sh) — kp{z) = —^^^6h. The bulk equations for the stroma 
keep their previous forms. 

Stress balance at the apical surface of the epithelium 
reads: 



i2rjid,v°)q6H 



0, 
0. 



(6) 
(7) 



The perturbation SH of the apical surface is determined 
by the boundary condition Vz\h+l+sh = loSH, which 
takes the form [kex-p{—H/l) — kQ]SH + Sv^ — lo5H to 
linear order. Stress balance at the epithelium-stroma in- 
terface reads: 

27]djvz - Spc = 2^idjuz - Sps + "fiq^Sh, 

= n{d:^Suz + djux). (8) 

Also at this interface, continuity of velocity and displace- 
ment yields {k—ko)Sh+Svz = ujSh and 5uz = 5h. Finally, 
the displacement vanishes at the bottom of the stroma: 
Sua\z=o = 0- The growth rate w is obtained by impos- 
ing the existence of a non-trivial solution to this set of 
linear equations. From this condition, wc obtain three 
relaxation modes for the system (Fig. [2]) . 

In the case where the stroma is treated as a viscous 
fluid, the previous equations need to be modified by re- 
placing the displacement Uq. by a velocity (w^) and the 
shear modulus /i by a viscosity (77s). In addition, the 
following boundary conditions are altered: the friction 
term in Eq. ([s]) and the condition Su^ = Sh at the 
epithelium-stroma interface are replaced by S,{Svx — Sv^) 
and Svl — cuSh, respectively. This results into two relax- 
ation modes (Fig. [s]). 

The number of modes that we get can be understood 
as follows. For the elastic stroma, the set of boundary 
conditions generates three modes because it contains the 
inverse relaxation rate uj three times: in the velocity- 
continuity conditions at both interfaces and in the tan- 
gential stress-balance equation at the epithelium-stroma 
interface. In the case of a fluid stroma, we loose the mode 
associated with the latter equality. 

It is instructive to look at the analytic expansions of 
these different modes in the limit of large wave numbers q. 
In this regime, the modes associated with respectively 
the epithelium-stroma interface and the apical surface 
decouple, since their characteristic decay lengths are of 
the order q^^, which is much smaller than H. For an 
elastic stroma, their expansions to constant order read: 
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(9) 



Among these expressions, only the one related to ojf 
can be positive, indicative of an unstable mode. It re- 
sults at the epithelium-stroma interface from a balance 
of the stabilizing surface tension and stroma resistance 
to deformations on the one hand, and the overall posi- 
tive cell-production on the other hand. This expression 
gives a necessary condition for the existence of an unsta- 
ble regime (ry(fc — fco) > /i). The condition = also 
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FIG. 3: Similar plots as those presented in Fig.[2]and with the 
same conventions, but in the case of a viscous stroma. The 
same default parameters are used, except for 77s = 10 kPa-s 
(instead of /x), 71 = 1 mN-m~^ and k = 0.9 d~^. 



FIG. 2: Relaxation modes w as a function of the wavenum- 
ber q for an elastic stroma, (a) The three relaxation modes 
are plotted for the following parameters: ry = 10 IVIPa-s, 
fj, = 100 Pa, 7i = 10 mN-m~\ 7a = 1 mN-m~^ (all estimated 
from [13]), k = 8.6 d"^ (see e.g. [U), ^ = 10^° Pa-s-m"^ 
(estimated from [T^), H — L — 300 fim and I = 200 fim (es- 
timated from [1]). (b) to (h) In each panel, the most unstable 
mode is investigated while one parameter is varied as com- 
pared with panel (a). The varied parameter is indicated at 
the top of each panel, and its different values directly on each 
graph. Plots are coded both in color as well as line styles. 



yields a leading-order expression for the upper crossover 
wavenumber from the unstable to the stable regime, pro- 
vided that this crossover occurs in the large-q domain. 
The expression for results from a balance of surface 
tension and cell production at the apical surface, and the 
one for uif from a balance of tangential stress and surface 
friction at the epithelium-stroma interface. Both expres- 
sions correspond to modes that are always stable in their 
region of validity. 

In the case of a viscous stroma, the potentially unstable 
mode reads: 
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(10) 



The second mode has an identical expansion to that of 
the elastic case, and the third mode is lost. 



Similar expansions can be obtained in the small-g 
regime, but the expressions to next-to-leading order are 
complicated and mix the different physical origins de- 
scribed above. In the case of an elastic stroma, two of 
the three relaxation rates diverge to minus infinity, in- 
dicative of the elastic resistance of the stroma to a uni- 
form compression. To leading order, they read: 
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The third mode however has a finite small-g limit, which 
reads fcexp {—H/l) — fcg. We can retrieve this expression 
by integrating the continuity equation &i q — Q over the 
height of the perturbed epithelium and to leading order 
in the perturbations. In the case of a fluid stroma, one of 
the modes has the same finite limit, which is consistent 
with the argument presented above. However, the other 
relaxation rate approaches zero as rather than infinity: 



[3g7a + 2£(7a + 7i)] 
677s 



(13) 



Therefore, as the system is also always stable at suffi- 
ciently small q, the relaxation time diverges in this case. 
This is because the relaxation here is associated with 
lubrication-like viscous flows over large distances in the 
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x-direction rather than elastic relaxation over short dis- 
tances in the ^-direction. 

These results show that the instability always occurs 
at finite wave vector. In Figs.[2]and[3| we analyze the be- 
havior of the most unstable mode as a function of the pa- 
rameters. We see that the interface is destabilized when 
either the epithelium viscosity rj, the cell-division rate 
k, or the thickness of the dividing region I is increased, 
because of a higher resulting shear stress [TH]. This is 
also true for the thickness L of the stroma in the elastic 
case, since a thicker stroma resists less to a given defor- 
mation. Increasing the other parameters has a stabilizing 
effect. This is intuitive for the elastic shear modulus of 
the stroma /z in the elastic case and the stroma viscosity 
77s in the viscous case, as well as for the surface tension 71 
in both cases. The parameter 7a in both cases as well as 
L in the fluid case have little influence on the dispersion 
curves (not shown for the fluid case). 

For a viscoelastic material with relaxation time t, we 
do not expect anything qualitatively different from the 
fluid or elastic cases to occur at large and intermediate 
wave vectors. In the small-g regime, as the relaxation 
rate goes toward a finite negative value in the case of 
an elastic stroma, it vanishes when the latter is fluid. 
Getting the correct behavior in the generic viscoelastic 
case would require a complete study. As a general fact, 
we expect the curves presented in Fig. [2] (resp. Fig. [s]) to 
be valid when ujt ^ 1 (resp. lot <^ 1). 

In this work, we have shown the existence of a hy- 
drodynamic instability of an interface between a viscous 
fluid with production terms and a viscoelastic material. 
The instability stems from the generation of viscous shear 
stress in the fluid due to material production. As such, 
this mechanism constitutes a new hydrodynamic insta- 
bility that has not yet been described. We propose that 
this effect provides a potential mechanism for the un- 
dulations at epithelium-stroma interfaces in vivo. Our 
analysis might explain why such undulations are more 
pronounced in neoplastic tissues [T]. Indeed, tumorous 
epithelial cells divide faster than healthy cells and in a 
larger domain away from the basement membrane I14j. 
The large-deformation regime of the instability might 
correspond to such fingering phenomena. It is commonly 
accepted that cancerous invasion requires the produc- 
tion of proteases that can degrade the basement mem- 
brane and remodel the extracellular matrix "14:]. Such a 
digestion could decrease the interfacial tension between 
the tissues as well as the elastic modulus of the stroma, 
thereby triggering the present instability. The digestion 
of the extracellular matrix is thus not an alternative to 
the mechanism proposed here, but one of its determi- 
nants. While proteases enhance the instability and al- 
low the growth of protrusions to proceed deeper into the 
stroma, we expect the physical forces driving this process 



to originate from the mechanism presented here. 

The undulation instability investigated in this work is 
potentially relevant for many biological systems in which 
interfaces of growing cell populations are present. For ex- 
ample, at interfaces between many tumors and healthy 
tissues, similar effects are observed [1]. More generally, 
we expect this type of instability to occur in all suffi- 
ciently viscous fluids with source terms. It would there- 
fore be interesting to conceive other systems that show 
the same type of instability, while being easier to char- 
acterize experimentally than living tissues. 

We thank A. Callan- Jones, M. Lenz and X. Sastre- 
Garau for many useful discussions. 
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